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1. Foreword and plumose of research. 


Determination in the field of the nature of ground surface material is 
frequently time and effort consuming, and indeed is not always technically 
possible. Thus the possibility of effecting such analysis by remote 
sensint ,y satellite or aeroplane) seem of great interest. 

At least theoretically, it may be possible to deduce the nature of the 
material forming a terrain from its thermal characteristics , by measuring, 
by remote sensing, for example, its surface temperature at two different 
and significant moments of the day, when temperature value is close to 
minimum and to maximum. In fact it can be shown [*l] , [Vj ; [Y]^ [uj 
that the value of total daily temperature excursion of a soil (maximum 
value minus minimum value) is connected to the product ^ . C . A 
which is characteristic of the material being examined, where is density, 
c is specific heat, \ thermal conductivity. 

On the whole this connexion is effectively recognisable only if the 
surface is composed of integral rocky material (i.e., has no surface 
stratum of degraded, fragmented or transported material), with a negligible 
humidity content, and not covered by any kind of vegetation ['LOG- 
As a hypothesis a rocky surface of the above kind was envisaged, consisting in 
blocks of three different materials (see par. 4, below) and an analysis 
was made of the possibility of identifying the dividing lines between such 
blocks and, where possible, to form a first approximate idea of the kind 
of material forming the surface of the terrain. 

This attempt at identification was carried out by simulation of the thermal 
dynamics of the soil surface during a day in which atmospheric conditions 
were average for the latitude in question (W. Europe) . 

As will be described in detail in par, 4, simulation was used to compute 
the values presented by the rocky surfaces consisting of blocks of three 
different materials, reference being made to spring or autumn months, at 
a latitude of about 4o°-50°. 
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The line of separation between zones of different materials is delineated 
by the range of temperature variation (difference between maximum and 
minimum temperature daily values', 


Numerical computation of t empera ture transi tor ies in a solid of finite 
dimension composed of vario us m ateria ls , using the Crank -Nicholson finite- 
"differences methods. 


f . I General , 

Widely-known analytical methods permit the solution of the differential 

equations governing diffusion of heat in solids 0], M, W, [9], [10], [uj 

These methods have certainly a theoretical interest since they make it 

possible to describe the physical phenomena with a general and organic approach 

nevertheless in practice they can only be used for the solution of rather 

simple problems, or as a quick way of obtaining rough preliminary results, 

since the weightiness of the calculations rapidly becomes prohibitive as 

soon as the boundary condition^ of the shape of the solid becomes complicated. 

The advent of computers produce a reduction in the importance of analytical 

methods for the solution of the differential equations, which have gradually 

. , , , mo died s , 

v'een replaced ov nusnur u,ni wiueh are much more easily implemented by 
the computer. A theoretical case of computation of temperature transitories 
under periodic regime was resolved using a finite-differences mathematical 
mo lei with time discretization according to the Orank-Nieholson method 
and solution of the system for computation of the unknown quantities by 
the Gauss -To. id el or surrolaxation method). For an exhaustive analysis of 
such finite-differences methods, see jp}j jidj |\lj • 


Numerical Case Study. 


A study was made of the temperature variations that would take place over 
a period of time on an island, rectangular in shape, and consisting of 
three different materials (fig. V, 

• c ic. 1 - riant of I Aland 
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It was assumed, as shewn in Tig. 1, that the separation surfaces between 
different materials were vertical planes having y*constant. 

The island itt supposed to be surrounded by the sea, the temperature of 
which, by a' first approximation, is taken to be at a constant value of 
d9 c;0 K throughout the day. In order to operate with the finite-differences 
method, the island was divided, on the x£ plane, into the network 
shown in Pig. 1.1, 


FIG. 1,1 - Dimension and subdivision in blocks of island 


Along the £ (depth) axis, temperature values are computed only for the 
first six metres, it being assumed that at such depths there is practically 
no diurnal temperature variation (it may be calculated with the numerical 
finite-difference? method?, or even taking an approximatenvalue from the 
analytical method, j~l] / that the daily temperature variation at such depth 
is of the order of 10 * °K for rocky materials. The temperature at 'o su 
depth is therefore assumed to be constant. 

On the vertical plane the solid is therefore divided into ID slices, 
each 0.? m thick (fig. 1.0), 


?IO. l.D - Subdivision of parallelepiped (island) 

blocks . 


into horizontal plane 


In fig, 1.1 is can be seen that the general (not borderline) mesh has 
the dimensions ^shown ill Big. I.2.I.J, 10 m. and. 20 m. j 


, while borderline meshes 

have a perpendicular' dimension at the every border of hQ era. (a borderline 
mesh is one belonging to the perimeter of the Island) . 


_ i. r ■ » i - 
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This differentiation va3 introduced to wake it possible to consider 
borderline meshes as being, at any given instant, in thermal equilibrium 
with the temperature of the sen ("waterline tone"). For such a border- 
line zone the order of magnitude of the perpendicular to border dimension 
was taken to be such that it would be permanently wetted by waves 
breaking on the coast, 

Boundary and initial conditions. 


The following boundary and initial conditions were assigned: 
Boundary conditions : 

VI 

a) Temperature constant, in time, at depth of t? m 


The temperature values at 6 m depth (values maintained constant in time) 
ore computed taking heat flow in a vertical direction to be constant 

and of value equal to that of the mean geothermal flow of the Earth, 

• • , u cal 

i.e., approximately l.s ± * 


cm 4 s 




It is further assumed that the initial surface temperature is uniform 
all over the island, the value being 09b°K. 

(This temperature value was chosen because it was observed, by the si- 
mulation of the transitory, to be the mean daily value reached by the 
soil surface temperature, with the values which were assigned in the 
example to the parameters belonging to the heat equation, when the 
initial temperature transitory in the simulation is fulfilled and the 
periodic stationary regime is reached) , 

Temperature at depth bOO cm. can easily be found from equation: 

A 

A, 

where 

K. 


£» * £u„ + 


.600 


cm 


( 1 ) 


I = heat conductivity of the i-th material (cal, 'em. s . °K) 
= surface temperature of material 
t^-yy* temperature at rOO cm depth 


<P 


= geothermal flow (cal/cm“ 
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Geothermal flow was token to be constant all over the island both 
because this' waa simpler, and in order to avoid the surface-temperature 
values computed during the transitory being influenced by the varia- 
bility of the geothermal flow. In fact, in the case- under consideration, 
the aim is to compute the influence of the various materials of which 
the island is composed, on the surface temperatures. 

b) Heat flow imposed nn surface of island . 

With the exception of the narrow "waterline" strip running along the 
perimeter of the island and 1*0 cm wide, in which, a3 already stated, 
surface temperature is token to be at the constant value of 293 °K (which 
is a reasonable supposition since the thermal inertia of water is much 
greater than that of the materials forming the island) , the net heat 
flow crossing unit ground surface in one second was supposed to be (fig. 2): 


FI! . 2 - Heat exchanges at the surface 


J (t) heat flow liy radiation from sun (eal/em w .s) in the surface 

An a first approximation for the sake of simplicity the heat flow is 
supposed to be sinusoidal, having a mean value of J . Hence the analytical 
expression of the heat flow reaching the surface by radiation from the 
sun will be of the following fonn: 

J (t) = J (l + sin (wfe) 

where eu i s pulsation of oscillation of solar radiation (l/s) 
w - 2 ~ ?T ’ where T » 2*» h. period 

T 


s 


. 5 . 
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value of radiation of Earth's surface towards space as given 
<hy the Stefan^Boltsmann lav (as a first approximation the 
temperature of the uni-verse is taken to he 0° K) , 


A 


7z' 


|.u»0 


geothermal flow crossing ground surface each second. 


<0 Temperatures imposed in points belonging to lateral surface (excluding 
upper and lower surfaces) of the paralielipiped (island). 


PIS. 3 - Section of island by plane x = const. (or y=const.) 


Fig. 3 represents a section of the island with a plane x 
y a constant. 


constant , or 


Fig. 1.1 gives a plan of the island. 

* lp " H ^*i» 2L» jk) denote the set of points bel-nging to the parallelepiped 
being studied (Fig, 3). 

fhe set of coordinates of the points forming the border of the island, 
to cm wide, shall have coordinates (x,y) belonging to a set S(x,y). 

Mow consider the points belonging to the lateral surface SL (x, y, s) 
of the parallelepiped and having Z < 6 m; i.e. SL (x,y,s) consists of 
the points with (x,y) 6 S (x,y) and with z 46 m, 
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. b , 


Such pointo belong to the set denoted by SL in Fig. 3.1 


I 


* ^ 

FIS. 3.1 - Section of island by planc^const. (or y= const.) 


Since there are reasons for believing that the areas belonging to SL 

which are closest to the surface (e.g. those with o / c d2m) ore affected 

• as regards K w 

by the influence of the sea temperaturh{snean temperature, it is supposed 
that the temperatures of such cones follow a pattern as a function of z 
of the type (fig. U). 


FI 3. 4 - Initial temperature values of parallelepiped (island) 


In Fig. 4 

initial temperature for cones belonging to SL. 



initial temperatures for cones composed of some material as above, 
belonging to S but not to 3L. 


t I . 
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It is supposed that at depths over 2a the influence of the temperature 
of the sea is' not felt, this influence tending to lower mean temperatures 
in zones immediately below the sea. 

To sum up, it was supposed that zones belonging to SL presented two 
temperature gradients in z\ one for 0 <z 42 m , and another for z > 2m (fig. 4). 
Lastly it was supposed (fixed potential boundary conditions), that the 
temperature of zones belonging to the set SL remained constant in time. 


Initial condition . 

Initial surface and deep temperature values. 

The initial value of surface temperature ^ is fixed. 


n3T(* x,o; = 2 s6 "K „/,„)£§■ 

\ * i 


% (x v J ■ -203 V ,, (*,/} 6 5 


(non-bunder zone) 
(border zone) 


for points having z > 0, <|) being geothermal flow and \ ^ thermal 
conductivity of the i - th material. We find 


(3) 

with Z/'O _ 

if point P (x,y,z) £ SL,or if point P belongs to SL and has z^2m. 

On the contrary, if point P (x,y,z) has 0<s 4 2wand belongs to SL, 
given that z^= 2m, we have the following initial temperature conditions 



%fo,z)-X(*,y l o) 


. 2 


(3.1) 


f 


. 8 » 


UNtVEHSITA OKOU #TUDI Dl MU, AN 
I8TITUTO 01 OKOFIIICA 


Thus (3) and (3*1) given an analytical expression of the conditions 
illustrated graphically in fig. k. 


3. Decryption of the finite-differences mathematical model used to 
simulate the temperature transitory of the island . 

On the (x-y) plane the island vas subdivided into a network of (8x1*0 
panes (fig) 

In direction z analysis was carried out for the first b m in depth, divi- 
ding this depth into I* cones (fig. 1.2). 

In total the parallelepiped was divided into 8 x I*< x 12 * I3M elementary 
parallelepipeds . 

Each elementary parallelepiped not belonging to the border set SI» thus has 
the dimensions (x, y, z respectively) of 20 tn x 10 m x 0.5 m (fig. 5). 


F1J. 5 - Dimensions of elementary parallelepiped (not 
belonjinj the border) 


the island as a whole having the dimensions ICO. 8 m x 120.8 rsx x 6 m. 

As regards the form of discretisation adopted, the following observations 
should be made: 

l) initially, since it was not easy to estimate a priori the magnitude 
of heat flows, simulation was carried out fo> a parallelepiped (island) 
consisting of elementary parallelepipeds of the same dimensions as the 
above but consisting, in directions x, y, s, respectively, of 12 x 20 x 20 
elementary parallelipiped (a total of 1 800) . 


I 
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It was however found that temperature variations in direction x,y 
being very "slow", and the temperature uci 5.5 m depth ’not varying 
significantly over time intervals of the order ot one week, it was 
possible to reduce the number of elements to 1 3^, as described above, 
without loss of precision. 

I) The computation programme was prepared in such a way as to be adaptable, 
by very 3mall end rapid modifications, to -un . with any number of 
elementary parallelipipeds having any given dimensions. The model is 
also suitable for use with a set of elementary parallelipipeds which at 
the limit might be composed each of different materials. 

i 

In the model, as a preliminary approximation, it was felt possible to 

neglect 3ome elements which undoubtedly affect surface temperature, the 

most important being (2] {V] : soil humidity (taken to be zero); evapora- 

' » * 

tion from soil; possible presence of degraded or fragmented material, or 
of surface vegetation; possible wind near ground surface; atmospheric 
temperature above the ground surface and air humidity; roughness of ground 
surface. 


Equations governing heat conduction in a solid (3~dimension model) . 

& 

The most, general differential equation governing heat conduction in a solid 
(in the absence of phase transformations) is [l] , [2]: 



( 4 ) 


(0 


J 


, 10 . 
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.2. feeing interior heat generation per unit volume and time 

i 

Having ...» apply CM to compute temperatures of a solid of finite 
dimensions, CM can be approximated with the corresponding finite- 
-differences equation [ii] . 

The reasoning whereby such an equation is obtained is as follows: 
considering the general case of a body fomed by elementary parallelipipeds 
having dimensions x, y, s differing from each other . 

As a generalisation, further supposing that the dimensions of each 
elementary parollelipiped differ from those of the others (fig. 6 ). 


. FI3. & - Blocks formed by two elementary parallelepipeds 
each vith different dimensions 


i 

V.'ith reference to the elements of fig. 6 compute the heat flow from 
cone N to zone C. 

In differential form the heat-flow equation is given by 



taking the joining line NC as x axis (fig. 6 ). 

The minus sign in (?) appears because heat flow versus is from the higher- 
-temperature to the lower-temperature sone. 

Defining as 
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R thermal resistivity (inverse of thermal conductivity^ 

The equivalent* thermal resistivity R NC between point N and point C 
is (in analogy to the sum of two resistences in series in an electrical 
circuit) . 

0 _ Eh Ay»+ R c Ay e (7) 

A>N + Ayc 

where R * resistivity of zone N, R_ that of zone C 
N ^ 

A.y^ A./c dimensions of parallelipiped in direction jr 
The equivalente thermal conductivity between zone N and zone C is (7)i 


v i - _Ay* t A y k - 

AnC ~ Rhc Rh^+Rc^c" 


A y H +Ay c _ A t1 A d Ay^ A)c 
-AA + Ay d t A c A/n 

Am /No 


Note that, since we are dealing with resistivity, (7) and (7.1) are 
obtained for two parallelipipeds of unit base area, hence for simplicity 


with 



FIG. 7 - Blocks formed by parallelepipeds of unit base area 


However, for the general case of fig. 7, we find that the total resistance 
between point N and point C has the value 



JL 

Ax c -tAz c 


Ah Ay c +Ac Ay^ 

2 Ah 


4. 

Ax c Az c 


(5.1) 


the total heat conductivity between point N and point C is, from (5.1 


(C _ ± . ^ . Ax.'Az 

T " _ R t „ A// c +\Ay N 



(5.2) 
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FIG. 3 - General blocks (c) with surrounding blocks 
(general dimensions) 
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Tbe boat passing between zone N and none C per unit time is, acccording 
to (frt and (5*2) (electrical analogues): 



(U)» : — Ax c Az t .^-^ ; 

i ■/ XmAy c +\a 



where A* c and Az^are the dimensions^ (with respect to x and a.) of the 
face perpendicular to y axis. 


In general, if we consider an elementary parallelipiped, C, surrounded 
by cones ( parallelepipeds) of dimensions in general differing from each 
other, the total heat flow c^, entering in C, each second from the 

- (fig S)) 


aur: 


ounding tones is, according to (5.2) 


q = j Q h A c AKt^a , (!fc-!fc) + hhi , 
' c [ \ h Ay c + Ac Ay„ X/ y f >My* ' y XrAxfi-XtA^ 



Aye A it $*-&)+ 

Ax c Ax», Az t |X C Az^j 



I n ae heat accumulating in one second inside the elementary parallelipiped G 
causer. , in time t , a variation of the temperature of G which is given 
by ('m 



(9) 


I 
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To the finite differences 


it 


is expressed as 


At 


For (8) and ($3 ve may thu3 write , taking coordinate to express the 
discrete value of x, i, that of y, and k that of a: 


~ fa ^ ± I JL ( 2 A-Ail 




At J c Uy 

Cu _ 

WtU 


+ 


X Ajk. ■‘"Aw Ay 
d. / 2 A A h 


2 


Ay 

i 




Ax ( ' AAx.+ A*Ax 


AX-t ^ L “"^ 1 0 


. jdSLti 


\A^,+Xu,' 

XXjVi 

9 AA^,+X; 4 X> 

AX 


251 + 2 


vw»i 




_____ , y 

t sr'xA^+v^ 7 2 z 


(40) 

In (10) in order to simplify notation the index in the expressions of Ax, Ay / 

Az Xd d 2" has been taken as understood 

when equal to i., j., k , 


Equation (10) makes it possible to compute the temperatures of the single 
• elementary parallelepipeds once the boundary conditions and the initial 
conditions of the system are known. The unknown temperatures t+&t) 

are computed from the known temperatures tX ( i ; j , k, t) at time t , 
using a system of equations of type (io). This system is resolved by the 

1 jll] 

The boundary heat-flow condition imposed at the surface of the solid, 
which is expressed in the differential form (2) , is expressed in the 

finite-differences form as follows. 


iterative surrelaxation method £d] 
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Let us consider a generic parallelipiped on the surface of the solid 
(s=0), hut not on the borderline, with geometrical center C - C \j( y t A? j 
(fig» 9) / ^ 


FIG. 9 ~ Heat exchanges in ^general surface blocks 


YOKiHi r. 

'*s g* wrv 


the heat equation for parallelipipeds belonging to the surface is expressed 
as : ^ 

A.AyAzffc A^l ) - JMt)* JW. Axdy- ) Ax-Av 

7 U At 1/ at 2. 1 t J 

T=t+ T 


.L^l Ax Ay+ Xm — • AyAz + L 3 AxAz 

l l 2 . 

'Y time variable (u) 

referred to point C / geometric center. 

Li discrete form (ll) is expressed by the following finite-differences 
equation, which is general, valid even in the case of z and \ 
varying w:»n j «pth ; 

d. d. 9 A A. 


& , i r wj» i 


AAz^+^Az Az 



\ + 1_ . 

/2U;-! /. 

^'v(r+At^"^/A 


z 

/ Ax 

X Ax 1-1 + ' 

2 

2- / 


z A 

\Axi+i AfA 


Z 2 


ZAA,'., 




\ A . i \ Av 


J) A/ UAy.^Ay 

. lX^y + X~ L b 
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In (12) tc simplify notation, the index is taken to he understood in the 
expression of when equal (respectively) to 

3 i 1 

A.* «,> A * 

It should be observed that the finite-differences method thus makes it 
possible to attain a reasonably simple solution of the heat equation 
also for a non-linear problem such as that studied here, furthermore in 

a region formed of elements of different materials and dimensions. 

~ ■ * 

The solution of such a problem by classical analytical methods would 
present monstrous difficulties of computing. Indeed, it may be stated 
that, since it is a non-linear problem, it cannot in practice be resolved 
by analytical traditional methods . 


H, Description of simulation of the thermal transitory .and analysis 
of results obtained. 


Figs. 1.1 and 1.2 show the region in study which consists of three diffe- 
ring materials. The width in the (x,y) plane of each of the three zones 
is respectively Ho ,8 m; Ho m, Ho, 8 m. 

The first zone consists of dolomite, the second of granite, the third 
of dry clay (i.e. with practically-negligible water content). 

It was supposed that the net solar heat flow J had a mean value of 15 m. 
.eai/cm 1 * s, which is a reasonable value for locations around Latitude Ho 0 

W ■ 

The following values were assigned to the parameters characteristic of the 
heat equations, They are roughly the mean for the materials in question 

[*U°] ■ 

Dolomite 

heat conductivity A = 10 cal/cm.s. *K 

*5 

density $ = 2,6 g/csir 

specific heat c - 0.23 cal/gk 
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Granite 

conductivity \ s 

7 . I0~ 3 cal/cm.s *K. 

density * 

2.52 g/cm^ 

specific heat c= 

0.20 cal/g.'K 

Clay 

conductivity * 

3,6.IO~ J cal/cm.s. *K. 

density * 

J.8 g/cm^ 

specific heat c® 

0.20 cal/g.'K 

The emissivity of 

the rock surface was 


all three materials ; furthermore it was supposed that of the heat irradiated 
toward the sky from the surface of the island, 25 # would be retained by 
the atmosphere, whilst 75# would pass beyond the barrier of the atmosphere 
and be irradiated toward the universe. 

The transitory was simulated starting from the initial condition described 
in par, 2, simulation being protracted for 110 h. It was found that after 
3 days (72h) the initial transitory could practically be taken to be 
concluded, a permanent periodic regime having been reached. 


From simulation it was deduced that heat exchanges in the x, y direction 
are very small compared with heat exchanges taking place at the surface 
due to heat flows in direction s_. 

In practice, for a given material, temperature differences between points 
of the surface reach at most 7/I00°K (in dolomite) ^Tab, a ^ Tab. a) at 
any given instant. 

tab a gives surface temperatures computed at a time approximately corre- 
sponding to a minimum of the mean value of the surface temperatures (<?6’h 


(°°) NOTE; the value of emissivity was taken to be equal for all three 
materials considered, since these three ‘types of rock effectively have rather 
similar mean emissivity values [is^ , though these may v^ry in 

particular cases over a certain range. 

It was thus felt that it would have been arbitrary to introduce any 
differentiation based on emissivity values for the three rocks under 
consideration . 
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after start of transitory). From top to bottom, zones are: dolomite, 
granite, clay.' 


TA3‘. A - Surface temperatures at time approximately 

corresponding to a maximum of mean surface temperalure 


Tab, a presents temperatures of the zones at a time approximatively 
corresponding to a maximum of the Bean value of surface temperature 
(IOT h after start of transitory). 

If we consider the central with respect to direction surface zones, 
i.e. close to the line W of fig. 10, 


, /V > -> 


Jurfaee 

spending 


temperature at time approximately corre- 
to a maximum of mean surface temperature 


we find the following temperature values, from top to bottom (i.e. from N 
to S) excluding border zones (TABLE l) , 


fAJ . 1 - temperature values along w line of fig. io 
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Cn the oth«r hand, if we consider the horizontal sections traced 
approximately ‘at the centre with respect to the £ axis of the single 

i 

hcmopreneous zones, in the E-W direction (lines 00, O' O', 0"0 M of fig. 
10. 1) still excluding border zones, we obtain the temperature values 
of TABLE II 


FI 3. 10 - Vertical axis (VV) of the Island 

t 


TABLE II 

TA3LS II - temperature values along OG, O'G', C"0" lines of 
fig. 10. 1 
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FI.'., 10.1 - Horizontal axis (00, Q'Q', 0"0") of each homogeneous 
zone of the Island 


As we ean be 3een in table I and II the passage from one material to 
another is marked by considerable variations of surface temperature values 
I and of 1, moving from North to South, while surface temperature and T 
variations are minimal the 3ome material both ir. the East-West and in the 
North-South directions (see also figs. 10 and lO.iJ* 


Thus in the case studied, each cone may for all practical purposes be 
considered to be nearly isolated' from the surrounding zones. 


FIG. -ti shows how, by thermal inertia measurements (i.e.,in practice by 
measuring the difference AT between the daily maximum and minimum tempera- 
tures), it is possible, in the simplified hypothetical situations considered 
• in this study, to identify materials of differing thermal inertias, It 
must be noted that in this study materials having widely-differing thermal 


inertias were chosen. 


In fact, for dolomite the value is £ - 0.0773 cal/ om“ 

cal /cm" \f5" *K; for clay O.OBbO cal/csT y s *K 
? = v/Xjc) • 


VT'K; for granite 
(thermal inertia 


/*• 


n * 
j* X 


Daily temperature variations 
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CONCLUSION, 


This report deals with the simulation of an ideal hut significant case of 
temperature transitory on the surface of an island consisting of different 
racks . 

Assuming certain simplifying hypotheses (listed in par, 3), numerical 

* 

simulation has shown the theoretical possibility of recognising, in the 
example dealt with, the transition from blocks consisting of one material 
to blocks consisting of another material having different thermal inertia 
by remote sensing of surface-temperature values, given the sensivity of 


present-day temperature sensors, 

Measurements t be made at least two different time;: of day; and if 
measurements are to be as significant as possible', ^Measurement should be 
made when the surface temperature values are close to minimum and the 
other when close to maximum. 

Xt has in fact been shown that the difference between maximum and minimum 
value of the temperature of ground surface during the day is linked to 
the thermal inertia value of the material of which the rock is formed. 
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